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Abstract 



We study the out-of-equilibrium evolution of an 0(2)-invariant scalar field in which a 
conserved charge is stored. We apply a loop expansion of the 2-particle irreducible effective 
action to 3-loop order. Equations of motion are derived which conserve both total charge and 
total energy yet allow for the effects of scattering whereby charge and energy can transfer 
between modes. Working in (l+l)-dimensions we solve the equations of motion numerically 
for a system knocked out of equilibrium by a sudden temperature quench. We examine the 
initial stages of the charge and energy redistribution. This provides a basis from which we 
can understand the formation of Bose- Einstein condensates from first principles. 

PACS numbers: 03.75.Nt, 11.30.Fs, 64.90.+b. 
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1 Introduction 

The formation of a Bose-Einstein condensate involves a significant proportion of the total charge 
in a system relocating to occupy the lowest energy state. This is achieved experimentally using 
dilute atomic gases where the conserved charge corresponds to the total number of atoms and the 
condensate is formed by reducing the temperature of the gas using laser and evaporative cooling. 
We can understand Bose-Einstein condensation (BEC) by considering the equilibrium state of a 



bosonic gas ^GIEIE]) but the question of precisely how the condensate arises as the gas responds 
to a change in its external conditions is a much more difficult problem (notable examples can be 
found in Refs.[SJ |SJ 0). This process is a very-many-body problem encompassing the properties 
of cold, trapped atoms and importantly, interatomic interactions. 

A theoretical understanding of condensate formation has recently been motivated by the pos- 
sibility of observing the spontaneous formation of defects 8 . The vortices already observed in 
rotating condensates could form spontaneously if the cooling process were to happen quickly 
enough. The mechanism by which this occurs has consequences not only for condensed matter 
physics but also for cosmology where it is predicted that similar defects may have formed sponta- 
neously in the early universe [5] . Experimental production of spontaneous defects in atomic gases 
would effectively allow us to test cosmological theories in the laboratory. 

In this article we shall consider the initial stages of the dynamical process by which charge 
and energy is redistributed between modes following a sudden drop in the ambient temperature. 
Recent progress has been made towards understanding such equilibration processes of quantum 
fields in Refs.[ini HH QH d d E]. In particular, the method of Refs.QI]! El El uses the loop 
expansion for the 2-particle irreducible (2PI) effective action where thermalization is observed at 
3-loop order |10| . Here we apply this method to a system containing a conserved finite charge. 

In the next section we outline our method for a relativistic scalar field theory invariant under 
0(2)-transformations (the non-relativistic theory emerges for low temperatures as is demonstrated 
in ^Hl)- This is one of the simplest models to have a conserved Noether charge, representing some 
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conserved quantum number. We derive equations of motion describing the evolution of the Green 
functions of the theory. In section |3| we consider the initial conditions needed to fully describe 
the system. For convenience we choose an initial equilibrium state. In section 0] we solve the 
equations of motion numerically. Due to computational constraints we work in one space and 
one time dimension. Though there is no phase transition in (l+l)-dimensions, trapped finite- 
size gases are known to exhibit interesting phenomena associated with a macroscopic occupation 
of the lowest energy state ^ E] ■ Here we examine the initial stages of the charge and energy 
redistribution which must form the basis by which BEC occurs. 



2 Dynamical equations 

To describe out-of-equilibrium behavior we use the Schwinger-Keldysh technique ^1 1191 I2()| . 
Given an initial distribution of states with density operator p at time U, the expectation of an 
operator O at some later time t is given in the Heisenberg picture by 

(6(t)) = Tr{p6(t)}. (1) 

In the path- integral formulation, the time contour C begins at time ij, evolves forward to time t, 
then returns to ti. 

We are interested in considering the response of a scalar field to a rapid drop in temperature. 
Assuming that the system is in equilibrium prior to the temperature quench, we may take the 
original equilibrium distribution as our initial distribution of states. Following the quench, as the 
system evolves with lower energy, it is forced to seek a new equilibrium. 

Since our initial state has non-zero charge, the grand canonical ensemble provides the obvious 
choice of density operator. Given a Hamiltonian H, charge Q, initial temperature T = 1/(3 and 
chemical potential p, this is given by 

/ 5 = exp{-/?(H-/iQ)}. (2) 

We are free to adjust the effective parameters of H in order to fix our initial distribution. The 
initial conditions are considered in more detail in section [3] 

The Hamiltonian density and charge operators for the self-interacting 0(2)-invariant scalar 
theory are 

Q = /"dx(7ri^ 2 -7T 2 0l), (3) 

where a = 1, 2 and repeated indices are summed over in the usual way. We proceed by writing 
out the Schwinger-Dyson equation for the 2-point Green function of this field (20j (we suppress 
field indices initially) 

(n x + m 2 )G2{xxi) + — G^xxxxi) — —i8(x — xi). (4) 
Diagrammatically, the 2-point Green function can be represented by a line as follows 

G 2 (xix 2 ) = xv — *x 2 (5) 

The 4-point Green function can then be written as a coupling expansion in terms of the 2-point 
Green function 

Xu ,X 2 X U ,X 3 X\, .#4 n Xi. .x 2 

G 4 (x 1X2X3X4) = + + - iX / dy \y +•••. (6) 

X3* »X4 X-f »X4 Xrf »X3 J X? • X4 
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Setting X2 = X2 = X3 — x we have 

Ga(xxxx\) = 3 — *x\ - iX J dy x^^y — 'X\ +■■■. (7) 

The first term on the right-hand side results from disconnected 2-point Green functions. Truncat- 
ing the series at this term constitutes a 'mean-field' linearization of Eq.QJ, known as the Hartree 
approximation. Since we wish to incorporate scattering of particles with the movement of charge 
and energy between modes, it is crucial that we keep at least the 0{\) term in the series expansion. 
It seems clear when we look at Eq.© that the last term is necessary if we wish to see scattering 
effects. 

Replacing field indices and neglecting the point-number label (all Green functions under con- 
sideration from now on are 2-point) we may rewrite the Schwinger-Dyson equation as 

(O x S ac + M^ c (x))G c b(x,y) + i dzY, ac (x, z)G cb {z,y) = -i6 ab Sc(x - y), (8) 

Jc 

where M is the effective mass and E is the non-local self-energy. Using the expansion of Eq. (J7J) 
and reinterpreting the algebraic form of the diagrams gives the effective mass as 

Ml b {x) = m 2 S ab + ^ [G cc {x, x)6 ab + 2G ab (x, x)} , (9) 
and the non-local self-energy as 
A x 2 



E a fc(x, y) = -2 l^-j [G cd {x, y)G cd (x, y)G ab {x, y) + 2G ac (x 1 y)G db {x, y)G dc {x, y)] . (10) 

We began with the Schwinger-Dyson equation and performed a coupling expansion of the 4-point 
Green function. The above result is equivalently obtained by starting with the 2PI effective action 
and truncating at 3- loop order as is shown in Refs. |llU13| (e.g. our Eq.© corresponds to Eq.(13) 
of Ref.|13)h In each of these references, a 1/N expansion is used in order to derive M and E. 
Although in principle our coupling expansion can be derived by truncating the I/TV expansion at 
C(A 2 ), this would require working to NNLO in order to encounter the (D(l/N 2 ) term on the right 
hand side of Eq. 1)10(1 . The weak-coupling limit of the NLO approximation is briefly considered in 
Ref.QS]. 

The problem we face is to solve Eq.© for some given set of initial conditions. The 2-point 
Green functions represent a time-contour ordered product of field operators 

G ab (x,y) = (Tc4> a (x)My))- (11) 

The usual procedure for solving the Schwinger-Dyson equation is to decompose the Green functions 
as [201 

G ab (x, y) = 6 c (xq - y )G a > b (x, y) + 6 c {yo - x Q )G< b (x, y) (12) 

such that 

G> b (x,y) = (Mx)Mv)), 

G< b (x,y) - (My)Mx))- (13) 
Now following the method outlined in Refs. |llM13| we take real and imaginary parts 

F ab {x lV ) = - [G> b (x,y) + G< b (x,y)] =ReG> b {x,y), 

p ab (x, y) = i [G> (x, y) - G<(x, y)] = -2ImG> b (x, y). (14) 
The function p is the spectral function and F is the symmetric propagator. 
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Similarly, we perform a decomposition of the self-energy 

Z ab (x, y) = d c {x - y )Y;> b (x, y) + 9 c (yo - x )T,< b (x, y), 
and again take real and imaginary parts 

Kb( x > v) = i K&fo y) - ^ab( x , y)] ■ 



(15) 



(16) 



We can now reexpress the Schwinger-Dyson equation in terms of F and p (cf. Eqs.(70) and 
(71) in Ref.[E]) 



(Px&ac + M 2 c )F cb {x,y) = - dz\ / dz Y, p ac {x,z)F cb {z,y) - / dz ^ c (x, z)p cb (z,y) 



(a x S ac + M ac )p cb (x,y) = - I dz<j / dz Y,P c (x,z)p cb (z,y) 

'Vo 



where now 



and 



M 2 b (x) = m 2 S ab + ^ [F cc (x,x)5 ab + 2F ab (x,x)} 
o 



(17) 
(18) 



yP 
^ab 



-21 i 



FcdFcd — -rPcdPcd ) F ab ~ ^FcdPcdPab 



/111 

"2 \F ac Fd b Fdc — -^FacPdbPdc — -^PacFdhPdc ~ ^PacPdbFdc 

FcdFcd — -^PcdPcdj Pab + 2F c dPcdF ab 



-2 ( PacFdbFdc + F ac pd b Fdc + F ac Fd b pdc ~ ^PacPdbPd< 



(19) 



Before attempting to solve these equations, we express them in momentum space. In general, 
given that J = J d d p/(2ir) d , we can write 



G(ti,f 2 ;x-y)= f e^-y^GihMw) 
Jp 



(20) 



due to spatial translation invariance. The self-energy terms T, F and S p take the form given 
above but factors in each term have momentum values pi, p 2 and (p — Pi — P2) respectively. 
Integrations are performed over pi and P2 and the self-energy carries a momentum label p. 
The dynamical equations become 

[(df 1 +p 2 )S ac + Ml(t 1 )]F cb (t 1 ,t 2 ;p) = - I f 1 dt'j:P ac (t 1 ,t';p)Fc b (t',t 2 ;p) 



dt'^ c {t 1 ,t';p)p cb (t',t 2 ;p) }, 



[(dl + p 2 )S ac + M 2 c (t 1 )]p cb {t 1 ,t 2 -p) = -\ dt'EP c {t u t';p)p eb {t',t 2 ;p)\, 



where the effective mass squared, expressed in the momentum basis is 

A " 



M 2 b [t) = m 2 5 ab 



G 



f F cc (t,t;p)6 ab + 2 f F ab (t,t;p) 



(21) 



(22) 
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It remains to solve these equations for some given set of initial conditions. Since our equations 
involve a weak-coupling expansion we should be careful to note that they are likely to give inac- 
curate answers should the coupling corrections become large. For this reason we should not stray 
too far from equilibrium. 

Finally, by taking expectations of the Hamiltonian and charge operators and substituting from 
Eq. Q for the 4-point Green function, we find the overall energy is given by 

E = [ E p (t) 
Jp 

| [d t d t 'F aa (t, t'; p)] t=t , =0 + (p 2 + m 2 + A j Fcc {t, t; p')) F aa (t, t; p) 
+ [zfij F ac {t,t; V 'Yj F ca (t,t;p) 

+ lJ o <tt' (ESc(t, P)F C a(t', t; p) - E£(f, f ; p) Pca (t', t; p)) | (23) 

and the overall charge is given by 

Q= I Q P (t) = I [d t F 12 (t,f;p) -d t F 21 (t,l/;p)] t=t , . (24) 
J p J p 

Each of these quantities can be shown both analytically and numerically to be conserved by the 
dynamical equations. However, charge and energy are indeed capable of exchange between modes. 




3 Initial Conditions 

The Schwinger-Dyson equations embody not only Heisenbergs equation of motion for the field but 
also the equal time commutation relations (ETCR) . It should come as no surprise that the ETCR 
are conserved by our dynamics so long as they are enforced by the initial conditions. The ETCR 
for our model are 

[0a(x,i),7r b (y,i)] = iS ab 6(x-y), 
[0 a (x, t), 6 (y, i)] = Mx, t), 7r b (y, t)) = 0. (25) 

These translate into statements about the spectral function 

[dtpab(t,t';p)] f=f , = 5ab, 
Pab (t, t- p) = [d t d vPah (t, t'; p)] t=t , = 0. (26) 

The ETCR are not sufficient to constrain the initial conditions. Other symmetry requirements 
which derive from the spectral function and symmetric propagator definitions Eq. l|14|l must also 
be satisfied 

F ab (t 1 ,t 2 ;p) = Fb a (t2,ti;p), 

Pab(ti,t 2 ;p) = -Pba(t2,h;p)- (27) 

(It can be checked that if these conditions are satisfied initially, then the dynamical equations 
preserve them for all later times.) Having satisfied the quantum commutation relations and other 
symmetry requirements, it remains to chose the initial conditions for these unspecified functions. 
Here we shall choose free field equilibrium values for some given temperature T — 1/(3, and 
chemical potential fi. Setting the initial time ti = we have 

F 12 (0,0;p) =F 21 (0,0;p) - 0, 



5 



F u (0,0;p)=F 22 (0,0;p) = — 

ZUJ 

[d t -F n (i,0;p)] t=o = [d t F 22 (t,0;p)] t=0 = 0, 
[d t F 12 (t,0;p)] t=0 = - [d t F 21 (t,0;p)] t=0 = - 



1 



1 



[dtd t >F ab (t,t';p)} t=t , =Q = LU 2 F ab (0,0;p), 



(28) 



where ui = \/p 2 + Mq , and Mq can be viewed as the effective mass squared prior to the temper- 
ature quench. 

Referring to Ea. (|24(l . the overall charge is given by 



Q = 



l 



l 



(29) 



This value is conserved throughout the subsequent evolution of Green functions. 



4 Numerical results 

In this section we present our numerical solutions to Eqs. (|21|l . 

Eqs. H21f) are integro-differential equations for F ab and p a b, each of which are functions of ti, 
t 2 and p. The equations involve coefficients which are themselves dependent on F a b(t,t';p) and 
Pab(t, t'; p) for t and t' up to and including t\ and t 2 , and for all values of p. The right hand side 
of Eqs l|21|) involve two integrations over internal momenta and one time integration. In addition 
these integrations must be performed for all external values of t\, t 2 and p. In light of this it 
would be computationally unfeasible to work in any greater than (l+l)-dimensions. Also, the 
range of the time integrations becomes greater as t\ and t 2 become greater and contrives to slow 
the calculations down at later times. This limits us to considering only early times. 

The (l+l)-dimensional continuum is replaced by a spatially periodic lattice. Integrations are 
performed by converting into discrete sums 



dp 1 



Na 



(30) 



N is the number of lattice points, a is the spatial lattice spacing, and Na = L is the length of the 
spatial dimension . Differential equations are solved by a basic leapfrog method. The time step 
h in units of 1/Minit is chosen to be 0.1 where Minit is the initial value of the effective mass 
following the quench. A corresponding choice of a = 0.4/Minit gives stability with h/a = 0.25. 
Choosing N = 100 we have a lattice size L — 40 in units of 1/Minit- With periodic boundary 
conditions the momenta are given by 



P -►-ra 



27T71 

lv~ 



I 



(31) 



Although these elementary techniques may provide a relatively poor approximation to the 
continuum, they are fast and provide stable results. Also, the conservation of energy, charge and 
other symmetries are transparent upon considering a single iteration of the dynamical process. 

Other parameters are chosen in units of Minit as T — 10, /i = 0.1 and Mq = 4. These 
correspond to a total energy of 172.92 and a total charge of 1.84. These values have no particular 
significance, they are chosen simply to represent the general behavior. Larger values of Mq provide 
a greater impulse, knocking the system further from equilibrium. Larger values of T reduce the 
impact of changing the mass parameter on the charge and energy distributions. The chemical 
potential must be less than Mq. For (3+l)-dimensional systems, the closer \i is to Mo, the closer 
the system is to the critical point. The coupling A = 0.5 is small in order to make sense of a 
coupling expansion yet large enough to see appreciable scattering effects relatively quickly. 



G 



M>(t) 




5 10 15 20 

Q„(t) 




5 10 15 20 t 

Figure 1: Evolution of effective mass squared, energy in the lowest mode, and charge in the lowest mode. 
The dashed line is the Hartree approximation, the solid line is the current approximation. 

FigQ] shows plots with the effective mass squared M 2 (t) fEa. (|22[l ). the energy in the lowest 
mode Eo{t) (see Ea.ip5)l). and the charge in the lowest mode Qo(t) (see Ea. lfSijl V The dashed 
line is the Hartree or mean-field approximation for comparison. In the Hartree approximation 
only terms of order A are kept in the dynamical equations. Solid lines represent solutions to 
Eqs.J2U- An important feature to observe is that the charge in the lowest energy state is constant 
in the Hartree approximation but is free to change at next order in A where scattering effects are 
incorporated. Oscillations in M 2 {t) and Eo(t) are damped out more effectively in the higher order 
calculation. 

FigHlshows the evolution of the symmetric propagator with time. For Fn(t, t; 0)(= F 2 2(t, t; 0)) 
we observe a much more effective damping of oscillations than in the Hartree approximation. We 
also observe an oscillatory decrease in correlations between the field at time t and the field at time 
0. This indicates that the system is moving towards some final outcome which is independent of the 
initial conditions and dependent only on the conserved quantities which characterize the system. 
It should be noted that the equations of motion may be reversed in time and initial conditions can 
always be recovered. A fixed point equilibrium solution in the future can in principle be approached 
arbitrarily closely but never reached JT]. There is zero equal time correlation between fields <pi 
and 4>2- 

FigUI shows the distribution of charge in the different momentum modes (Eq.JSJJ)). The upper 
figure indicates the initial charge distribution (see Ea.l|29jl). The lower momentum modes store a 
greater amount of charge. As the system evolves we see that charge moves out-of the region p > 1 
and into the region p < 1. The change in charge at each momentum value is shown at times 0, 10, 
and 20. This process fits our intuition where we would expect that charge moves into the lower 
energy states as the temperature decreases. 

The energy distribution of Eq. is given in Fig0J The upper figure gives the initial energy 
distribution. In the continuum we would expect the energy to increase quadratically with p and 
without bound. Our lattice approximation flattens out the energy spectrum at the upper limit 
of momentum whilst the spectrum is approximately quadratic at low momenta. Following the 
quench we see a movement of the energy from the region p < 1 into the region p > 1 . The change 
in the energy distribution becomes less oscillatory as time increases. This is a reflection of the 
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Figure 2: Evolution of the symmetric propagator. The dashed line is the Hartree approximation, the 
solid line is the current approximation. 



fact that oscillations in both the effective mass and the equal time propagator are damped out as 
time increases. 

FigEl displays the increase in spatial correlations as a result of the drop in temperature. The 
correlation function 

F u (t,t;x) = f e ,;px Fii(M;p) (32) 
Jp 

shows a growth in amplitude and a growth in the size of the region in which the field is correlated. 
These regions can be interpreted as domains in which the field has strong correlations |21j . The rate 
at which these domains grow has consequences for the formation of defects in higher dimensional 
systems undergoing a symmetry breaking phase transition |22) . In our simulations it is found that 
the correlation function quickly grows to a new stable function which then remains approximately 
constant. 

In addition, in our numerics we observe conservation of total charge, conservation of total 
energy, and conservation of commutation relations. 



5 Conclusions 

We have considered the quantum field dynamics of an 0(2)-invariant scalar field theory with a 
conserved charge. We have argued that in deriving our equations of motion, it is crucial that we 
work to at least 3-loop order in the 2-particle irreducible effective action. This is the lowest order at 
which the effects of scattering are included, allowing for the movement of charge between modes. In 
lower order calculations, the charge in each mode is artificially confined to that mode. Since BEC 
is essentially characterized by the charge distribution, in order to consider the formation process, 
the charge in the system must be able to redistribute in response to some external stimulus. 
More generally, to understand any equilibration process involving finite conserved charge we must 
account for the transfer of charge between modes. 

We have considered a (l+l)-dimensional system. Although this does not allow us to model a 
condensation process, we are able to examine the charge and energy distributions as they evolve 
in time. For a sudden drop in temperature modelled by a sudden change in the effective mass 
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-0.01 

Figure 3: Charge distributions at a series of times. 

squared, we observe an expected general movement of charge from higher to lower momentum 
modes. This is accompanied by a movement of energy from lower to higher momentum modes. 

A (3+l)-dimensional version of this calculation remains a desirable goal for the future. An 
improvement of the algorithm in order to run the calculation to much later times would also be 
welcome. 
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